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Abstract. For linear and fully non-linear diffusion equations of Bellman- 
Isaacs type, we introduce a class of monotone approximation schemes rely- 
ing on monotone interpolation. As opposed to classical numerical methods, 
these schemes converge for degenerate diffusion equations having general non- 
diagonal dominant coefficient matrices. Such schemes have to have a wide 
stencil in general. Besides providing a unifying framework for several known 
first order accurate schemes, our class of schemes also includes more efficient 
versions, and a new second order scheme that converges only for essentially 
monotone solutions. The methods are easy to implement and analyze, and 
they are more efficient than some other known schemes. We prove stability 
and convergence of the schemes in the general case, and provide error estimates 
in the convex case which are robust in the sense that they apply to degenerate 
equations and non-smooth solutions. The methods are extensively tested. 



1. Introduction 

In this paper we introduce and analyze a class of monotone approximation 
schemes for fully non-linear diffusion equations of Bellman-Isaacs type, 

(1.1) u t - inf sup \ L a > p [u](t, x) +c a 'P(t,x)u + f a 'P{t,x)\ =0 in Q T , 

(1.2) u{Q,x)=g(x) in R N , 
where Q T := (0,T] x and 

L a ^[u]{t, x) = tr[a a '"(t, x)D 2 u{t, x)] + b a ^{t, x)Du(t, x). 

The coefficients a a ^ = \u a ^a a ^ T , b a ^ ', c a <P, f a ^ and the initial data g take 
values respectively in S , the space of N x N symmetric matrices, ~R N , M, R, and 
R. We will only assume that a a ^ is positive semi-definite, thus the equation is 
allowed to degenerate and does not have any smooth solutions in general. Under 
suitable assumptions (see Section [5]), the initial value problem (|l.l[l - (|1.2[) has a 
unique, bounded, Holder continuous, viscosity solution u. This function is the 
upper or lower value of a stochastic differential game, or, if A or B is a singleton, 
the value function of a finite horizon, optimal stochastic control problem |25j . 
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We introduce a family of schemes that we call Semi-Lagrangian (SL) schemes. 
They are a type difference-interpolation schemes and arise as time-discretizations 
of the following semi-discrete equation 

Ut- inf sup \Ll' [l h u](t,x) + c ' 13 {t, x)u + f a ^ (t, x)\ = in X ,l x(0,T), 

where L"'^ is a monotone difference approximation of L a ^ and Z/, is a monotone 
interpolation operator on the spatial grid X . For more details see Section [5J 
Typically these scheme are first order wide-stencil schemes, and if the matrix a a ^ 
is bad enough, the stencil has to keep increasing as the grid is refined to have 
convergence. They include as special cases schemes from [16l HH [12], more 
efficient versions of these schemes, and a new second order compact stencil scheme. 
There are two main advantages of these schemes: (i) they are easy to understand 
and implement, and more importantly, (ii) they are consistent and monotone for 
every positive semi-definite diffusion matrix a a ^ = ^a a '^a a '^ T . The last point is 
important because monotone methods are known to converge to the correct solution 
[3J, while non-monotone methods need not converge [21] or can even converge to 
false solutions [2"3"] . 

Classical finite difference approximations (FDMs) of (|1.1[) are not monotone (of 
positive type) unless the matrix a a ^ satisfies additional assumptions like e. g. being 
diagonally dominant [18]. More general assumptions are given in e.g. [3 [14] but 
at the cost of increased stencil length. In fact, Dong and Krylov [14] proved that 
no fixed stencil FDM can approximate equations with a second derivative term 
involving a general positive semi-definite matrix function a a '@ . Note that this type 
of result has been known for a long time, see e.g. [20[ I12j . Some very simple 
examples of such "bad" matrices are given by 

(i 4 H^in[0,l] 2 , K "f) for .4 = 3 =[0,1], I - 
\±x 1 x 2 x{ ) \a(3 (3 Z J L J ' \Du\ 2 

and these type of matrices appear in Finance, Stochastic Control Theory, and Mean 
Curvature Motion. The third example leads to quasi-linear equations and will not 
be considered here, we refer instead to [12"] . 

To obtain convergent or monotone methods for problems involving non-diagonally 
dominant matrices, we know of two strategies: (i) The classical method of rotating 
the coordinate system locally to obtain diagonally dominant matrices s a ^, see e. g. 
Section 5.4 in |18j . and (ii) the use of wide stencil methods. The two solutions 
seem to be somewhat related, but at least the defining ideas and implementation 
are different. Both ways lead to methods that have reduced order compared to 
standard schemes for diagonally dominant problems. But it seems to us that it is 
much more difficult to implement the first strategy. 

In addition to the wide-stencil methods mentioned above, we also mention the 
method of Bonnans-Zidani [5] which is not an SL type scheme. Schemes for other 
types of equations related to our SL schemes have been studied by Crandall-Lions 
[T2] for the Mean Curvature Motion equation, by Oberman [35] for Monge- Ampere 
equations, and by Camilli and the second author for non-local Bellman equations [8]. 
The terminology SL schemes is already used for schemes for transport equations, 
conservation laws, and first order Hamilton-Jacobi equations. In the Hamilton- 
Jacobi setting, these schemes go back to the 1983 paper [9] of Capuzzo-Dolcetta. 
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The rest of this paper is organized as follows. In the next Section we explain 
the notation and state a well-posedness and regularity result for (|l.lj) - (ll.2| . The 
SL schemes are motivated and defined in an abstract setting in Section [31 and in 
Section|4]we prove that they are consistent, monotone, L°°-stable, and convergent. 
We provide several examples of SL schemes in Section [5j including the linear in- 
terpolation SL (LISL) scheme. This is the basic example of this paper, and it is a 
first order scheme that can be defined on unstructured grids. 

Our SL schemes make use of monotone interpolation, and higher order interpola- 
tion is not monotone in general. But for essentially monotone solutions, we can use 
monotone cubic Hermite interpolation (see Fritsch and Carlson [17] and Eisenstat, 
Jackson and Lewis |15j ) to obtain new second order schemes called monotone cubic 
interpolation SL (MCSL) schemes. These new schemes are defined in Section [6j 
and in contrast to the LISL schemes they are compact stencil schemes. Note well 
that in the special case of first order HJ-equations with monotone solutions, these 
schemes are consistent, monotone, second order schemes! To our knowledge, this is 
the first example of a monotone scheme which is more than first order accurate in 
the HJ-setting. 

We discuss various issues concerning the SL schemes in Section [7] We compare 
the LISL scheme to the scheme of Bonnans-Zidani [5 and find that the LISL scheme 
is much easier to understand and implement, and in general it is much faster. 
However, on bounded domains the LISL scheme will in some cases over-step the 
boundaries and some ad hoc solution has to be found. This problem is avoided by 
the Bonnans-Zidani scheme at the cost of being less accurate near the boundary 
than in the interior. Finally, we explain that the SL schemes can be interpreted 
as collocation methods for derivative free equations, or as dynamic programming 
equations of discrete stochastic differential games or optimal control problems. 

In Sections [H [9] and Appendix [Bj we produce robust error estimates for convex 
equations (i. e. B is a singleton in (jl.lj) ). These estimates are obtained through the 
regularization method of Krylov and apply to degenerate equations, non-smooth 
solutions, and both the LISL and MCSL schemes. Finally, in Section[Tni our meth- 
ods are extensively tested. In particular we find that the LISL and MCSL schemes 
yield much faster methods to solve the finance problem of Bruder, Bokanowski, 
Maroso, and Zidani [6|. 



2. Notation and well-posedness 

In this section we introduce notation and assumptions, and give a well-posedness 
and comparison result for the initial value problem - (|1.2|l . 

We denote by < the component by component ordering in R M and the ordering 
in the sense of positive semi-definite matrices in S w . The symbols A and V denote 
the minimum respectively the maximum. By | • | we mean the Euclidean vector 
norm in any M. p type space (including the spaces of matrices and tensors) . Hence 
if X e R N ' P , then \X\ 2 = £\ ,. \X tJ | 2 = tr(XJf T ) where X T is the transpose of X. 

If w is a bounded function from some set Q' C Q x into either K, R M , or the 
space of N x P matrices, we set 

\w\ = sup \w(t,y)\. 
(t,y)eQ' 
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Furthermore, for 8 £ (0, 1], we set 

r i \ w (ti x ) ~ w ( s i U)\ ii 'i r i 

Ms = sup j i i, njks and \w\s = Mo + Ms- 

Let C b (Q') and C°' 5 {Q'), 5 £ (0, 1], denote respectively the space of bounded con- 
tinuous functions on Q' and the subset of Cb(Q') in which the norm | • \$ is finite. 
Note in particular the choices Q' — Qt and Q' = M. N . In the following we always 
suppress the domain Q' when writing norms. 

For simplicity, we will use the following assumptions on the data of (|l.l|l - (|1.2|l : 
(Al) For any a £ A and f3 £ B, a a ^ = ±a a < f3 <r a < f3T for some N x P matrix cr Q A 
Moreover, there is a constant K independent of a, (3 such that 

\g\x + W a '?\i + \b a <P\i + Ic^K + l/^K < K. 

These assumptions are standard and ensure comparison and well-posedness of 
(|1.1|) - (|1.2|) in the class of bounded x-Lipschitz functions. 

Proposition 2.1. Assume that (Al) holds. Then there exists a unique solution u 
of (|1.1[) - (|1.2[) and a constant C only depending on T and K from (Al) such that 

Mi < c. 

Furthermore, if u\ andu-i are sub- and supersolutions of satisfying ui(0, •) < 
1*2(0, ■), then i*i < i*2. 

The proof is standard. Assumption (Al) can be relaxed in many ways, e. g. 
using weighted norms, Holder or uniform continuity, etc. But in doing so, solutions 
can become unbounded and less smooth, and the analysis becomes harder and less 
transparent. Therefore we will not consider such extensions in this paper. 

By solutions in this paper we always mean viscosity solutions, see e.g. [TT1 [25] . 

3. Definition of SL schemes 

In this section we propose a class of monotone approximation schemes for (|l.ip - 
(|1.2[) which we call Semi-Lagrangian or SL schemes. This class includes (parabolic 
versions of) the "control schemes" of Menaldi [19] and Camilli and Falcone [7] and 
the monotone schemes of Crandall and Lions [12]. It also includes SL schemes 
for first order Bellman equations [9l [16] , and it allows for more effective versions 
of these schemes as discussed in Section For a motivation for the name, we 
refer to Remark 17.21 For the time discretization we propose a generalized mid- 
point rule that includes explicit, implicit, and a second order Crank-Nicolson like 
approximation. Note that the equation is non-smooth in general, so the usual way 
of defining a Crank-Nicolson scheme would only give a first order scheme in time. 

The schemes will be defined on a possibly unstructured family of grids {Ga*,Ax} 
with 

G = GAt.Ax — {(tn, Xi)}neN ,ieN = {t n }nefio X ^Ax, 

for At, Ax > 0. Here = to < t\ <••■<*„ < t n+ i satisfy 

max At n < At where At n = t n — i„_i, 

n 

and Xax = {xi}ieN is the set of vertices or nodes for a non-degenerate polyhedral 
subdivision T Ax = {Tf x } je ® of R N . For some p £ (0, 1) the polyhedrons T 3 = Tf- X 
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satisfy 

int(Tj n Ti) = 0, II Tj = R N , pAx < sup{diam B Tj } < supjdiamT.,} < Ax, 

where diam is the diameter of the set and is the greatest ball contained in Tj. 

To motivate the numerical schemes, we write a = (eri, a%, a m , erp) where 
a m is the m-th column of a and observe that for k > and smooth functions 6, 



I 1 P 

-tT[aa T D 2 c/>(x)] = - ]T tr[a m alD 2 4>(x)} 



m=l 
p 

El 4>(x + ka m ) - 24>{x) + 4>(x - ka m ) | n[ ,2\ 
2 k 2 j ' 

m— 1 

& ^(q ; ) = ^ + fc ^ ) "^ ) +0(fe 2 ) 

_ 1 <t>{x + fc 2 b) - 2<j>{x) + 0(s + fc 2 6) 2 

These approximations are monotone (of positive type) and the errors are bounded 
by jgP\a\Q\D 4 (f>\ k 2 and i 1 1«| § | Z? 2 0| fc 2 respectively. To relate these approximations 
to a grid G, we replace <fi by its interpolant I<p on that grid and obtain 

^ r Tr>2i/ m 1 ( 1 9){x + ka m )-2(2^)(x) + (2(j))(x-ka m ) 

m— 1 

1 + fc 2 6) - 2(J0)(x) + (T^Kg + k 2 b) 

If the interpolation is monotone (positive) then the full discretization is still mono- 
tone and represents a typical example of the discretizations we consider below. 

To construct the general scheme, we generalize the above construction. We now 
consider general finite difference approximations of the differential operator L a ^[(/>] 
in (jl.ip defined as 

, qi , rQ/W , ^ </>(t, x + y«f + (t, x)) - 2<j>{t, x) + </>(t, x + yZf-fr x)) 

(3.1) L k ' p [(f>}{t,x) — , 

i=i 

for k > and some M > 1, where for all smooth functions </>, 

(3.2) \L a /y>] - L^[<f>}\ < C(\D<j>\ + ■■■ + \D^\o)k 2 . 

This approximation and interpolation yield a semi-discrete approximation of (jl.ip , 

U t - inf sup \L?^[lU](t,x) +c a - p (t,x)U + f a -P(t,x)X = in (0,T)xI Al , 
and the final scheme can then be found after discretizing in time using a parameter 

0e[o,i], 

(3.3) 8 At JJ? = inf sup \ [1U 6 ^ 1+6 + c ^-^uf' n + f^-i+e^ 
mG, where U? = U(t n , Xi ), f^~ l + e = f^{t n ^+9At n , Xi ), ... for (t n , Xl ) S G, 

W(t,») = ^ 3?) ~^~ At>ai) , and #.» = (i_^-i 
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As initial conditions we take 

(3.4) Uf=g(xi) in X Ax . 

Remark 3.1. For the choices 8 = 0,1, and 1/2 the time discretization corresponds 
to respectively explicit Euler, implicit Euler and midpoint rule. For 9 = 1/2, the 
full scheme can be seen as generalized Crank-Nicolson type discretizations. 



4. Analysis of SL schemes 



In this section we prove that the SL scheme p.3|) is consistent and monotone, 
and we present L°°-stability, existence, uniqueness, and convergence results for 
these schemes. In Section [8] we also give error estimates when B is a singleton and 
equation (|1.1[) is convex. 

For the approximation and interpolation X we will assume that 

M 



(Yl) 



(II) 
(12) 



i=i 

M 



EKf + 

i=l 
M 

E 

i=i 

M 

E 



Vki 



a, 3, 

Vki 



3 ,.a,P,+ 



4 

j=l£/fc,i 



5 3 a,/3, 



.4 ce,/3, 



= 0(k% 
= 0{k 4 ). 



There are X > 0, r G N such that \(l(f>) - <p\ < K\D p (j>\ Ax p for 

all N 3 p < r and all smooth functions 4>. 

There is a non-negative basis of functions {wj(x)}j such that 

(I(j))(x) — ^^(f>(xj)wj(x), Wi(xj) = Sij, and Wj(x) > for all i,j G 



Under assumption (|Y1[) . a Taylor expansion shows that L^' 13 is a second order 
consistent approximation satisfying (|3.2p . If we assume also (|TTj) , it then follows 
that Lfr[T<j>] is a consistent approximation of L a,/3 [0] if ^§ > 0. Indeed 

where |££ ,/3 '[</>] - L"'^]! is estimated in (|3T2|) . and by ^TTJ> , 



|L^ET^]-L^[0]|<C|ir^o 



Remark 4.1. Assumption (|Y1|) is similar to the local consistency conditions used in 
[18] . The 0(k A ) terms insure that the method is second order accurate as k — > 0. 
Convergence will still be achieved if we relax 0(k i ) to o(/c 2 ) as k — > 0. 



Remark 4.2. An interpolation satisfying (}I2j) is said to be positive or monotone and 
preserves monotonicity of the data. Note that such an interpolation X</> does not 
use (exact) derivatives to reconstruct the function <f>. Typically I will be constant, 
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linear, or multi-linear interpolation (i.e. r < 2 in (jTTj) ) since higher order interpo- 
lation is not monotone in general. For later use we note that from (|TTj) and fl2|) it 
follows that 



(4.1) r > 1 =»- y^^Wj(x) = 1 and r > 2 V] XiWj( 



x) = x. 



Now we prove that the scheme <|3.3|l — (|3.4p is consistent and monotone. The 
scheme is said to be monotone if it can be written as 

(4.2) supinf {Bff '"'"I/] 1 - V Bjf'^UJ 1 - V Bff •"•""Vf- 1 - F/'M = 

in G, where £}? ! >' 3 > n > m > o and F" ,,3 ' n does not depend on {J. 

Lemma 4.1. 4s SM me (|TT|) 7 dHJ, and (|YT|) /io/d. 

(a,) T/ie consistency error of the scheme (|3.3p is bounded by 

l 1_2 ^l|J. I \ -t i f~1 1 A+2/1, I , I I , I ni I , I n 2 



\cj> tt \ At + Cl At 2 (|<fe| + 1</> ttt\o + \D(f> tt \ + \D 2 <p 

tt\o) 

Ax r 



+\D r <f>\o=r + (\D<t>\o + ■■■ + \DU\o)k 2 ^ 

(b) The scheme (|3.3p is monotone if the following CFL condition holds, 
-M 



(4.3) {1-9) At 



a,P,n-l+8 

k 2 ° l 



< 1 and 9Atc"'^' n 1+8 < 1 for all a, f3,n,i. 



Remark 4.3. By parabolic regularity D 2 ~ <9 t which means that e.g. \D 2 (j> tt \ is 
proportional to | <^> ttt | o - When 9 — 1/2 ("Crank-Nicolson"), the scheme (|3.3[) is 
second order accurate in time. 

Proof. It is immediate that the scheme (13.31) is consistent with (jl.ip with a trun- 
cation error bounded by 

Jl_^|0 tt | o Ai-4|0 ttt |oAi 2 + sup {|£^[^]-i^Cr^"]| } 

* «> a,/3,n 



sup 



for smooth 0. By p} and |L Q ^[/<"] - L^' [l(f> e ' n ]\ can be bounded by 

C|£> r 0| o -^- + C(|^|o + • ■ ■ + |^V|o)fc 2 , 
while \L a >P[4> n - l+e - 4> 8 ' n }\ is bounded by 

At 2 9{l - 0)sup |L Qj3 [(Mo < CAt 2 9(l - 0){\D<t>tt\o + \D 2 fr t \ \. 

a,f3 1 J 

Finally, | c «A™-i+0(^n-i+0 _ < C9(l - 9)At 2 \(j) tt \^. Hence part (a) follows. 

To prove part (b), we note that since J2i Wi = 1 we have 

L^mt, ■)}(t n - 1+0 ,x j ) = lff' n ~ 1+e [<t>(t, Xi ) ~ 4>{t, Xj )] , 
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where 



a,P,n-l+8 
hi 



M 

£ 

1=1 



Wi(xj + y^f' + (t n - 1+g ,Xj)) + Wi(xj + y%'f' (t n -\ +e ,Xj)) 



This quantity is non-negative by (12), and JDi^i 1+ " = p- by (|H 



M Jm |J2| 5 and 

J2i u>i = l. Therefore the only non-zero coefficients in (|4.2|) at (t n ,Xj) are 



B?f' n ' n = f + 9At n (# - j?'".**-^* _ c «,/J,»-i+«) 



i,i 

.a,/3,n,n — 1 
id 



6AtJ 



a,/3,n-l+9 



,a,p,n-l+$ 



a,/3,n-l+t 



), 



where j ^i. These coefficients are positive if (|4.3p holds. 



□ 



Existence, uniqueness, stability, and convergence results are given below. In the 
following, we denote by c a '@' + the positive part of c"'' 3 . 

Theorem 4.2. Assume (Al), (HQ), I©, fYT]) , and (|4~3|) . 

(oj There exists a unique bounded solution Uu of (|3.3p — (|3.4[) . 



fftj TTie solution U h of (|3~3|) - (|3~4)) is L°° -stable when 26 At sup a 



< e 2sup ° 



• + |o«. 



|Sf|o + t n sup|/ 



a,/3] 



a,/3 



(c,) iT^ converges uniformly to the solution u of (|l.l|) - (jl .2p as At, fc 



A.i 



< 1; 



0. 



Proof. Existence and uniqueness of bounded solutions follow by induction. Let 
i = t„ and assume U n ~ 1 is a known bounded function. For e > we define the 
operator T by 



TEH 1 



*7J 



e • (left hand side of (02])) for all jeZ 



,1/ 



Note that the fixed point equation U n = TU n is equivalent to equation (|3.3[) . By 
the definition and sign of the ^-coefficients we see that 



Tin 1 - Tin 1 



{ [l - £ (1 + At n 9( 



< sup ■ 

< (1 -e[l- A* n 0sup|c' 

a,(3 



„ce,P,n-l+t 



)) M 



■sAtJ^\U n -u: 



A,/3,+ | 



if e is so small that l-£(l+At6>(S-c' 



a,/3,n-l+t 



)) > 0and£(l-At n 6»sup ai/; 



lo)< 



1 for all j,n,a,(3. Taking the supremum over all j and interchanging the role of 
U and U proves that T is a contraction on the Banach space of bounded functions 
on Xai under the sup-norm. Existence and uniqueness then follows from the fixed 
point theorem (for U n ) and for all of U by induction since U° = g is bounded. 



A similar argument using (|4.2[) proves L°° - stability for 9Atsup a * 



|C/ r ' 



< 



a,/3,+| 



1 + (l-6>) At supple 1 
f-0Atsup Qi/3 |c a A+| o 



< e 2su P Q ^ I 



lot™ 



|.g|o + t„ sup|/ 



\u n -\ 



o,0\ 



At n SUp|/ 



a,P\ 



<*,P 
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In view of this estimate, convergence of Uh to the solution u of (|1.1|) - (|1.2|) follows 
from the Barles-Souganidis result in [3]. □ 

5. Examples of SL schemes 

5.1. Examples of approximations L^' 13 . We present several examples of approx- 
imations of the term L a,/3 [cj>] of the form [<j>], including previous approximations 
that have appeared in pTH]. \19\ [7l [12] plus more computational efficient variants. 

1. The approximation of Falcone [16] (see also [9]), 

b *,P D(j) ra % + 

corresponds to our \ik = Vh, y^'^ — k 2 b a -^ . 

2. The approximation of Crandall-Lions |13j . 

p , l(f>(x + ka a ^) - 210(x) +l(j){x - ka^ 13 ) 



itrK^ T ^]«^: 



2 rj ^ 2fc 2 

corresponds to our if 'j 3 '* = ±ka"' 13 and M = P. 

3. The corrected version of the approximation of Camilli- Falcone [7] (see also [19]') . 

itr[(T Q ' /3 (J a ' /3T D 2 (/ ) ] + b"'' 3 ^ 

,A J0(a: + Vhaf P + - 21<j ) {x) + lcf>(x - Vho-*' 13 + p a <P) 

2h ' 

corresponds to our L^ ,/3 if fc = vft, y^f ,=t = ±kaf' fi + ^fe"'' 3 and M = P. 

4. The new approximation obtained by combining approximations 1 and 2, 

itr[cr Q < / V' /3T L> 2 0] + b a ' P D<j> 

~ P + ^ 2P ' 

corresponds to our if y^'f^ = ika"'^ for j < P, y^ p^ = k 2 b a ^ and 
M = P + 1. 

5. The new, more efficient version of approximation 3, 

itr[<7 Q < / V' /3T P> 2 0] + b^Dcj) 

^ Z^(i + kaf j ) - 214>(x) +l(f>(x - kaj' ) 

~ ^ 2P 

i=i 

l<t>{x + ka a / + k 2 b a ^) - 21(f>(x) + l<j>(x - ka a / + k 2 b a ^) 
+ 2k 2 ' 

corresponds to our if y%j' — ±ka"'^ for j < P, y^p :± = ±fccp /3 + k 2 b a ' 13 
and M = P. 

Approximation 5 is always more efficient than 3 in the sense that it requires 
fewer arithmetic operations. The most efficient of approximations 3, 4, and 5, is 4 
when <j a 'P does not depend on a, but b a 'P does, and 5 in the other cases. 
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5.2. Linear interpolation SL scheme (LISL). To keep the scheme (|3.3|) mono- 
tone, linear or multi-linear interpolation is the most accurate interpolation one can 
use in general. In this typical case we call the full scheme (I3.3|) - (|3.4p the LISL 
scheme, and we will now summarize the results of Section [4] for this special case. 

Corollary 5.1. Assume that (Al) and (lYlj) hold. 

(a) The LISL scheme is monotone if the CFL conditions (I4.3[) hold. 

(b) The consistency error of the LISL scheme is 0(\1 - 26\At + At 2 + k 2 + ^-), 
and hence it is first order accurate when k — 0(Ax 1 ^ 2 ) and At — O(Ax) for 9 =/= i 
or At = OiAx 1 / 2 ) for 9 = \. 

(c) If 29 At sup Q ^ |c"'^ ,+ |o < 1 and (|4. 3[) hold, then there exists a unique bounded 
and L°° -stable solution Uh of the LISL scheme converging uniformly to the solution 
u of (01)- (Ti~2j) as At,k,^^ 0. 

From this result it follows that the scheme is at most first order accurate, has 
wide and increasing stencil and a good CFL condition. From the consistency error 
and the definition of L^ the stencil is wide since the scheme is consistent only if 
Ax/k — » as Ax — > and has stencil length proportional to 

max J^f-|V|^f + l k 

I := — t ~ > oo as Ax — ► 0. 

Aa; Aa; 

Here we have used that if (Al) holds and a ^ 0, then typically Uuf' ~ k. Also 

note that if k = Ax 1/2 , then I - Ax~ 1/2 . Finally, in the case 9 ^ 1 the CFL 
condition for (13. 3|) is At < Ck 2 ~ Aa; when k = 0{Ax 1 / 2 ), and it is much less 
restrictive than the usual parabolic CFL condition, At = (9(Aa; 2 ). 

6. A SECOND ORDER SL SCHEME FOR MONOTONE SOLUTIONS 

In this section we introduce new second order SL schemes of the form (|3.3|) - 
(|3.4|) for non-degenerate grids of tensor product type. These schemes are based 
on monotone cubic Hermite interpolation |17[ 1 1 5 j , and are consistent for monotone 
solutions of the scheme. We will call these schemes monotone cubic interpolation 
SL schemes or MCSL schemes in short. 

To define monotone cubic Hermite interpolation, we start by considering a ID 
function <f>. For each sub-interval [a;i,a;i+i], i € Z, we construct a cubic Hermite 
interpolant 

fulfilling 



{I(j>){x) = C + Ci(x - Xi) + c 2 (x - Xi) 2 + c 3 (x - Xi) 3 



(14>)(xi) = <pi, (14>)'(xi) = di, 

(lcj>)(x i+1 ) = 4> i+ i, (l(j>y(x l+1 ) = d i+ i, 

where <fn = 4>{xi) and di is an estimate of the derivative of 4> at Xi. It follows that 

(6.1a) c = (f>i, ci = di, 

,„ , v 3Ai - d l+ i - 2di 2Aj - d l+1 - d t 

(6.1b) c 2 = r , c 3 = 



Aa; ' J Aa; 2 



SL SCHEMES 



11 



where A; = ^'"^ — ■ To obtain a fourth order accurate interpolant, must be at 
least third order accurate. We will use the symmetric fourth order approximation 

(6.2) d l = _ , zeZ. 

However, the resulting interpolation is not monotone. Necessary and sufficient 
conditions for monotonicity were found by Fritsch and Carlson [T7] (see also |24| ) : 
If Aj = 0, then monotonicity follows if and only if di = dj+i = 0, and if 

di jo 
ai = - and ft = — , 

then monotonicity for Aj ^ follows if and only if (cuj, ft) eM = A4 e U Mb where 

M e = {(a, 13): (a - l) 2 + (a - l)(/3 - 1) + (/3 - l) 2 - 3(a + /? - 2) < 0}, 

M b = {(a, ft) : < a < 3, < < 3}. 

Eisenstat, Jackson and Lewis [T5] give an algorithm that modifies the derivative 
approximation di such that the above conditions are fulfilled, and for monotone 
data the resulting interpolant is a C 1 fourth order approximation. We will only 
consider C° interpolants, and in that case their algorithm simplifies to: 
Step 1 Compute the initial di using (|6.2[) . 

Step 2 Compute A,. If A, ^ compute on and ft, else set a, = ft = 1. 
Step 3 Set 014 := max{ai,0} and ft := max{/3j,0}. 
Step 4 If (a,, ft) ^ A^, modify (aj,ft) as follows: 

• If a, > 3 and ft > 3, set = ft = 3, 

• else if ft > 3 and a^ + ft > 4, decrease ft such that (aj,ft) € 9.M, 

• else if ft > 3 and + ft < 4, increase such that (a,, ft) € 

or a, = 4 — ft , in the last case subsequently decrease ft until 

( ai ,pi)edM, 

• else if ai > 3 and c^+ft > 4, decrease on such that (aj, ft) € dM, 

• else if aj > 3 and ai + ft < 4, increase ft such that (ai, ft) G <9.M 

or (3i — 4 — a j , in the last case subsequently decrease aj until 
KftJe&M. 

In each sub-interval [a;,, iCj+i], we replace di by ajAj and d,-+i by ftAj in (|6.ip . 
Multidimensional interpolation operators are obtained as tensor products of one- 
dimensional interpolation operators, i. e. by interpolating dimension by dimension. 

Lemma 6.1. The above monotone cubic interpolation is always monotone and 
satisfies (jl2|) . If the interpolated function is strictly monotone between grid points, 
then (llTj) holds with r = 4 and the method is fourth order accurate. 

Proof. Assumption (|I2|) holds by construction. The error estimate follows from 
[15] . since the above algorithm coincides with the two sweep algorithm given there 
when n = 1 interval is considered. In [15) it is proved that this algorithm gives third 
order accurate approximations to the exact derivatives and hence the cubic Hcrmitc 
polynomial constructed using this approximation is fourth order accurate. □ 

Remark 6.1. Carlson and Fritsch [lOj constructed an alternative monotone bicu- 
bic interpolation algorithm for R 2 . We are not aware of any work on high order 
monotone interpolation on unstructured grids. 

By the Lemma l6~Tl and the results in Section [4] we have the following result: 
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Corollary 6.2. Assume that (Al), (jYljl hold, and that for all h 6 (0, 1), i/iere 
exists a bounded solution Uh of the MCSL scheme such that TUh is strictly x- 
monotone between points in the x-grid X& x . 

(a) The MCSL scheme is monotone if the CFL condition (|4.3[) hold. 

(b) The consistency error of the MCSL scheme is 0(\1 - 20|Ai + At 2 + k 2 + 

and hence the scheme is second order accurate when k = O(Ax) and At — 0(Ax 2 ) 
for0^ \ or At = O(Ax) for 9 = §. 

(c) If26Atswp a p |c a, ^|o < 1, then the solution Uh is unique, L°° -stable, and con- 
verges uniformly to the solution u of (jl.ip - (jl.2[) as At, k, ^§ ► 0. 

Remark 6.2. The MCSL scheme is always monotone, but for non-monotone solu- 
tions the scheme is not consistent when k — O(Ax) since the consistency error then 
is 0{At + k + 4§-). Moreover, it is easy to see that the MCSL scheme has strictly 
monotone solutions (between grid points) whenever the collocation equation (|7.ip 
(see Section [7]) has strictly monotone solutions (between grid points). To prove 
such type of results one can use comparison principle arguments, and we refer to 
Appendix IA1 for results concerning equation Since (|7.ip satisfies the compar- 

ison principle (the proof is essentially a simplified version of the proof in Appendix 
IB[) . the argument proving Lemma I A . 1 1 shows that its solutions are monotone under 
assumption (A2) in Appendix lAl Under this assumption, existence of a monotone 
solution follows from Theorem 14.21 in Section 0J 

Remark 6.3. It is well known that consistent monotone methods for first order 
equations are at most first order accurate in general. However, the MCSL scheme 
is an example showing that second order consistent monotone schemes are pos- 
sible if the solutions have special structure: When cr a,/3 = 0, the MCSL scheme 
is a monotone second order scheme for a first order equation when solutions are 
monotone. Experts we have talked to seem to be surprised by this fact. 

7. Discussion 

7.1. Comparison with the scheme of Bonnans-Zidani (BZ). In the paper [5] 
(see also [4J [6] ) Bonnans and Zidani suggest an alternative approach to discretize 
non-linear degenerate diffusion equations. Their idea is to approximate the diffu- 
sion matrix a"'' 3 by a nicer matrix a^'^ which admits monotone finite difference 
approximations. For every k € N they find a stencil 

5 fc c{^Ki,..,We^: 0<max|&|<fc, i=l,...,N} 
and positive numbers such that 

This new diffusion matrix gives a diffusion term that can be decomposed into 
a linear combination of directional derivatives and these are again approximated 
by central difference approximations, 

tr[o^DV] «tr[ O ^U a 0] = £ a^D^ « £ <f A,0, 
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where D| = tr[££ T £> 2 ] = (f • D) 2 and 

A^w(x) = {w(a; + £Aa;) — 2ui(a;) + — £Aa;)}. 

This approximation is monotone by construction and respects the grid. In two 
space dimensions, a?'^ can be chosen such that 

| a a,/3_ a «,/3| = 0(k~ 2 ) (cf. [4]), and 
then it is easy to see that the consistency error is 

0{k~ 2 + k 2 Ax 2 ). 

When b 01 ' = 0, the BZ scheme can be obtained from (|3.3p by replacing our 
L?'^ by the above Bonnans-Zidani diffusion approximation. This scheme shares 
many properties with the LISL scheme, it is at most first order accurate (take 
k ~ Aa; -1 / 2 ), it has a similar wide and increasing stencil, and it has a similar good 
CFL condition At < Ck 2 Ax 2 (~ Ax when k ~ Aa; -1 / 2 ). To understand why the 
stencil is wide, simply note that k by definition is the stencil length and that the 
scheme is consistent only if k — > oo and kAx — > 0. The typical stencil length is 
k ~ Aa; -1 / 2 , just as it was for the LISL scheme. 

The main drawback of this method is that it is costly since we must compute the 
matrix for every x, t, a, (3 in the grid. The LISL scheme is easier to understand 
and implement and runs faster. Later we will see numerical indications that LISL 
runs at least 10 times faster than the BZ scheme on some test problems. The BZ 
scheme has the advantage that it is easy to modify to prevent it from leaving the 
domain (accuracy is then reduced or monotonicity is lost). It is less natural to do 
this for the LISL scheme. However, in many problems it is not necessary to do any 
modification near the boundary as we will see below. 

The MCSL scheme in the typical case when k = Ax, is a second order accurate 
and compact stencil scheme having the usual (not so good) CFL conditions for 
parabolic problems At ~ k 2 — Ax 2 . It is far more efficient than the other two 
schemes as can be seen in Section [lOl but it is only guaranteed to converge when 
the computed solutions are essentially monotone. The other two schemes "always" 
converge. We have summarized our findings in the table below. 





order 


CFL 


wide 
stencil 


boundary 
conditions 


convergence 


efficiency 


BZ 


1 


At ~ Aa; 


yes 


OK 


always 


worst 


LISL 


1 


At ~ Aa; 


yes 


special 
treatment 


always 




MCSL 


2 


At ~ Aa; 2 


no 


OK 


monotone 
solutions 


best 



Table 1: Comparison between the BZ, LISL, and MCSL schemes. 



7.2. Boundary conditions. When solving PDEs on bounded domains, the SL 
schemes (and the BZ schemes) may exceed the domain if they are not modified 
near the boundary. The reason is of course the wide stencil. This may or may 
not be a problem depending on the equation and the type of boundary condition: 
(i) For Dirichlet conditions the scheme needs to be modified near the boundary or 



14 



DEBRABANT AND JAKOBSEN 



boundary conditions must be extrapolated. This may result in a loss of accuracy 
or monotonicity near the boundary, (ii) Homogeneous Neumann conditions can be 
implemented exactly by extending in the normal direction the values of the solution 
on the boundary to the exterior, (iii) If the boundary has no regular points, no 
boundary conditions can be imposed. In this case the SL schemes may not leave the 
domain because the normal diffusion tends to zero fast enough when the boundary 
is approached. Typical examples are Black-Scholes type of equations. 



7.3. Interpretation as a collocation method. The scheme (|3.3[l - (|3.4[l can be 

interpreted as a collocation method for a derivative free equation, this is essentially 
the approach of Falcone et al. [16j E] ■ The idea is that if 

W Ax (Qt) — {u : u is a function on Qt satisfying u = lu in Qt} 

denotes the interpolant space associated to the interpolation I, equation (|3.3[) can 
be stated in the following equivalent way: Find U € W Ax (Qt) solving 

in i \ s tjtl • r f r cc.SrfrQ.mn— 1+9 , a.S.n— 1+9 fr&.n . ea.B.n— 1+6 \ • 

(7.1) S Atn Ul = mf 8 upiL k ' p [U ]i +c i ' p +f i ' p ' > m G. 

In general W Ax can be any space of approximations which is interpolating on the 
grid Iai! e. g. a space of splines. We do not consider this case here. 



7.4. Stochastic game/control interpretation. In general the scheme (|3.3|) - 
(|3.4[) can be interpreted as the dynamical programming equation of a discrete 
stochastic differential game. We will now try to explain this in the less techni- 
cal case when B is a singleton and the game simplifies to an optimal stochastic 
control problem. 

Assume that (Al) holds, and for simplicity, that c a (t,x) = and the other 
coefficients are independent of t. Then it is well-known (cf. [25 ) that the (viscosity) 
solution u of equation (|1.1|) - (|1.2|) is the value function of the stochastic control 
problem: 

(7.2) u(T-t,x)= min e\ [ f a{s) (X s ) ds + g(X T ) 

where A is a set of admissible ^-valued controls and the diffusion process X s = 
j£t,x,a(-) j g con strained to satisfy the SDE 

(7.3) X t = x and dX s = a a(s) (X s ) dW s + b a[s) ds for s > t. 

This result is a consequence of dynamical programming (DP) and (II. ip is called the 
DP equation of the control problem (|7.2p - (|7.3|) . Similarly, the schemes (I3.3[) — (I3.4[) 
are DP equations (at least in the explicit case) of suitably chosen discrete time and 
space control problems approximating (|7.2[) - (|7.3p . We refer to [JH] for more details. 

We take the slightly different approach explored in [21 HFJ HH H] to show the 
relation to control theory. The idea is to write the SL scheme in collocation form 
(|7.ip and show that (|7.ip is the DP equation of a discrete time continuous space 
optimal control problem. We illustrate this approach by deriving an explicit scheme 
involving L% as defined in part 4 Section ISTD Let {t — 0, ti, ...,tM — T} be 
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discrete times and consider the discrete time approximation of (17.21) - ()7.3|) given by 

M-l 



(7.4) u(T -t m ,x) = min E \ V f a <° (X k ) At k+1 + g(X M ) 

k=m 

(7.5) X m =x, X n = X„_i +cr Q "(l„_i)fc„C« + 6 a "(X„_i)fc2 77„, n > m, 



where fc„ = y/(P + l)At n , Am C A is an appropriate subset of piecewise constant 
controls, and = (£ ra ,i> ■ ■ ■ , £n,p ) T and r] n are mutually independent sequences of 
i. i. d. random variables satisfying 

i , ((fn,l,...,£n ) P,»k) = ±e i ) ^ 2(p 1 + ^ if j€{l,...,P}, 

P(j,€n,l, ■ ■ ■ , £n,P) ?7n) = e P+l) = p 

(ej denotes the j-th unit vector) and all other values of (£ n ,i, ■ ■ ■ , £n.p,?7n) have 
probability zero. Here we have used a weak Euler approximation of the SDE coupled 
with a quadrature approximation of the integral. Note that X rs X and u ~ u when 
At is small. By DP 



n-l 



u(T-< m ,.T)= min E\Y^ f ak {X k )At k+1 +u{T -t ni X n ) 



for all n > m, 



k—m 



and taking n = m + 1, s M -m = T - t m , As m = s m - s m _i, fe m = k M _ m , and 
evaluating the expectation using (|7.5p . we see that 

u(sM-m,^) = 

r ^ 2 1 

min I f a (x)As M - m + -^rrrT~ L t^_ , [u](sM-m-i,x) + 2(sAr- m -i,i) k 

I ST -\- X M m— 1 J 

where L£ is as in Section |5. II part 4. If we subtract u(sm— m— from both sides 

and divide by Asm-™ = 'p^' 1 , we find (|7.1[) with = 0. 

In [7], a similar argument is given in the stationary case for schemes involving 
the L£ of part 3 Section IQI In fact it is possible to identify all L^s appearing 
in Section [5. II with DP equations of suitably chosen discrete time continuous space 
control problems. However assumption (|Y1[) is not strong enough for this approach 
to work for the general defined in Sections [3] and [4j 

Remark 7.1. A DP approach naturally leads to explicit methods for time dependent 
PDEs. But implicit methods can also be derived using a trick. Discretize the PDE 
in time by backward Euler to find a (sequence of) stationary PDEs and use the 
DP approach on each stationary PDE. For stationary problems the DP equation is 
always implicit, so the result is an implicit iteration scheme. 

Remark 7.2. By the definition of and (lYljl . x + yf' k can be seen as a short time 
approximation of (|7.3|) . Hence the scheme (|3.3[) tracks particle paths approximately. 
In fact by the above discussion we might say that the scheme follows particles in 
the mean because of the expectation. For first order PDEs, schemes defined in this 
way are called SL schemes by e. g. Falcone. Moreover, in this case our schemes will 
coincide with the SL schemes of Falcone |16j in the explicit case. This explains why 
we choose to call these schemes SL schemes also in the general case. 
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8. Error estimates in the convex case 

We will derive error bounds in the case when B is a singleton and is convex. 
It is not known how to prove such results in the general case. Here and in the 
following, we do not indicate the trivial dependence any more. For simplicity we 
also take a uniform time-grid, letting G — At {0, 1, . . . , Nt} X X& x in this section. 
Let Qa*,t := At {0,1,..., Nt} X R w and consider the intermediate equation 

(8.1) S At V n (x) = 

inf \L%[V d > n }(t,x)+c a (t,x)V e > n (x)+f a (t,x)} in Q At , T , 

(8.2) V(0,x) = g(x) in R N . 

The first step is now to find a bound on \U — V\. 

Lemma 8.1. Assume that (II), (12) and the CFL condition (|4. 3|) hold and that 
sup„ \V n \i < C v . IfV solves l[8TT j) - ([8T2j) and U solves (pH]) - (jiT!) . then 

At 

\U-V\<C— in G. 

fx, 

Proof. Let W = U — V and subtract the equation for V from the one for U to find 

W? < + At sup f Ll[lW 6 ^- 1+e + cl n - 1+e W°' n 

ueA 



Assuming V n has p bounded derivatives, we rearrange the equation and use (|TT|) to 
see that 

(l + 0At(^-c-r))wp 
< W?- 1 + At sup |#(L£[ZVnr 1+e + T^W, n 

a£A ^ ^ " 

+ (1-9) {L^TW^ 1 ]^ + cf^-^Wr 1 ) } 

A T rAp 

+ 2AtKsup\D rAp V n \ — — in G 

n & 

with c^j +e Wl l = sup Q c l a '"" 1+ V™. By the CFL condition g3]), the coefficients 
of the above inequality are all non-negative. Hence since W n < iW^lo := sup 4 |Wj n |, 
we may replace W n by |W n |o on the right hand side. Moreover, since IlT^lo = 
|V^"|o and L^lH 7 ™^] = 0, the upper bound on the right hand side then reduces to 

M ~ Ar rAp 

(1 + At(l - 6)C c )\W n -% + 6At^\W n \ + AtK — p — , 

where C c = max a |c" ,+ |o. The same bound holds if wc replace W by —W, and 

hence we can conclude that 

M M ~ Ar rAp 

(1 + Ai0(— -C c ))\W n \ < (l + Atil-e^lW^lo + eAt—lW^o + AtK—^. 

Since W° = in X Ax , an iteration then reveals that 

m=0 
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(Y2) 



when At is small enough. Since V n is Lipschitz (p = 1), the lemma follows. □ 

Next we want to estimate \V — u\ when u solves (|l.ip - (|1.2[) . This can be done 
using the regularization method of Krylov if we can find suitable continuity and 
continuous dependence results for the scheme. These results rely on the following 
additional (covariance-type) assumptions: Whenever two sets of data cr, b and a, b 
are given, the corresponding approximations L^,y^'f and L^,^'^ 1 in (|3.ip satisfy 

M 

+ i&1 - m + & < 2k 2 (b« - 6°), 

i=l 
Af 

Er a,+ ^ i ct . — ^ a. — i . r~a.+ , ~q.— ^ ~q. — i 

bfe.i ® y k ,i + Vk,i ® ] + ® + «k,i ® «k,i ] 

-mm ®fiM ®i/k,< + +y k ; l ®Vkj\ 

< 2/c 2 (a Q - 0((j Q - ct q ) t + 2fc 4 (6 a - 6 a )(6 Q - 6 a ) T , 

when cr, 6, are evaluated at (t, x) and ct, b, y k are evaluated at (t, y) for all t, x, y. 
In Section [5] we will prove the following error estimate. 

Theorem 8.2. Assume that B is a singleton, that (Al), (lYTj) . l|Y2]) . and i/je CFL 
conditions (|4.3p ZioZd, and i/iai fc G (0, 1) and At < (2ko A 2/ci)" 1 . If u and V are 
bounded solutions of [fTT ]) -([P ]) and (|53 ]) -([83 ]l . i/ien 

<C(|l-26'|Ai 1 / 4 + At 1 / 3 + /fc 1 / 2 ) zn Q At , x . 

It also follows from the regularity results in Section [5] (see Proposition [HH]) that 
\V n \i < 2Ct, so by Lemma IQ1 and Theorem 18.21 we have the following result. 

Corollary 8.3 (Error Bound). Under (II), (12), and the assumptions of Theorem 
Iffl if u solves (fl~l"j) - (fT2"j) and U solves (|33|) - §tty , then 

\u-U\ < \u-V\ + \V -U\ < C(\ 1-20IA* 1 / 4 + At 1 / 3 + k 1 / 2 + ^-) in G. 

k 

This error bound applies to both the LISL and MCSL schemes, and it also holds 
for unstructured grids. If the solutions are more regular, it is possible to obtain 
better error estimates. But general and optimal results are not available. The 
best estimate in our case is 0(Ax 1 ^ 5 ) which is achieved when k = 0(Ax 2 ^ 5 ) and 
At = 0(k 2 ). Note that the CFL conditions (|3~3j) already imply that At = 0(k 2 ) 
if 9 < 1. Also note that the above bound does not show convergence when k is 
optimal for the LISL scheme (k = 0(Ax 1/2 )) or the MCSL scheme (Jfe = O(Ax)). 

Remark 8.1. These results are consistent with results for LISL type schemes for 
stationary Bellman equations. In fact if all coefficients are independent of time and 
c a (x) < — c < 0, then by combining the results of [7] and [T], exactly the same error 
estimate is obtained for the solution of a particular stationary LISL scheme and 
the unique stationary Lipschitz solution of . 



9. Proof of Theorem 18.21 

We start by an existence and uniqueness result. 

Lemma 9.1. Assume that (Al), (|Y1|) . and the CFL conditions (|4.3ll hold. Then 
there exists a unique solution Uh £ Cb(Qr,At) of (|8.1|) - (|8.2p . 
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The proof is similar to (but simpler than) the proof of Theorem 14.21 with the 
modification that the fixed point is achieved in the Banach space C(,(R ) instead 
of the space of bounded functions on Xj\ x . 

We will now give a result comparing subsolutions of (|8.ip to supersolutions of 

S At U n (x) = inf {L%[U e > n ](t,x) + c a (t,x)U^ n + f a (t,x)\ 

ctG-4 l J 

(9- 1 ) inRV>l, 
U(0,x) =g(x) in R N , 

where L£ is the operator defined in (|3.ip , (jYip , (|Y2I) when a a , b a are replaced by 
b a ~b a . 

Theorem 9.2. Assume that (Al), (|Y2) . (|4T3]) /io/d /or feot/i ([EI]) and ([9T]). // 

E/ G C(QT,At) is a bounded above subsolution of (|8.ip and [/ G C(Qx,At) & bounded 
below supersolution of (19. ip . i/ien /or a/Z G (0, 1), At < (ko A 1 J x,y £ WL N , 
nG {0,1,...,A T } 7 

U(t n ,x)-U(t n ,y) < R ko (t n )\(U(0,-)~U(0,-))+\ Q 

+ RkoiU^Rktitn)^ + t n L)\x - y\ 

+ t n sup [|(/ - /)+| + H fco (t n )(|C7| A |L>| )|c- c| ] 

+ ^ /2 2A' T sup [\b-b\ + \a-a\ ] 

where R k (t) = 1/(1 - fcAi) t/At ; A T < R ko (T)R kl (T)(L a + TV), 

u = isii v l = (|c a i! v \c a \i)(\u\ a if/| ) + irii v iri l5 

fco =sup|c Q '+|o, fci =2sup{| - Q |? + |6 Q |2 + l}. 

a a 

Remark 9.1. The function R k (nAt) = 1/(1 - ftAt)" satisfies 

5/\tRk{t n ) = kR k {t n ), 

R k (Q) = 1, and R k {t n ) < e 2kt " when At < ±. 

This is a key result in this paper, and the proof is given in Appendix |Bl In the 
stationary case, results of this type have been obtained in [UH] for simpler schemes. 
The result is a joint uniqueness result (take (er, b, c, /, g) — (a, b, c, /, g)), continuous 
dependence result (take x = y), boundedness, and rr-Lipschitz continuity result: 

Corollary 9.3. Under the assumptions of Theorem \ 9.2l if k G (0,1) and At < 
(2k A 2ki)~ 1 , then any bounded solution U G CbiQr.At) of (|8.ip satisfies 

ft) \U(t n ,-)\o <e 2kot "\g\ +t nS up a \r\ , 

(ii) \U(t n ,x) - U{t n ,y)\ < e 2(fco+fe!)t„( Lo + tnL) \ x _ y \ t 

where the constants, which are defined in Theorem \9. l A are independent ofk, At, Ax. 



Proof. Part (i) follows from Theorem [921 and Remark 19 . 1 1 since U = satisfies (|9.ip 
with (& a ,b a ,d a J a ,g a ) = (cr a , b a , c a , 0, 0). Part (ii) follows by taking U = U and 
x^y. □ 
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Now we extend the scheme (|8.1|) to the whole space Qt- One way to do this 
and to obtain continuous in time solutions, is to pose initial conditions on [0, At) 
by interpolating between g[x) and U(At,x) where U is the solution of (|8.1[) - (|8.2p . 

(9.2) S At V(t, x) = inf {it [V e (t, •)] if ,x) + c a {t e , x)V e (t, x) + f a , x) } 

in (At, T] x R N , 

(9.3) V(t,x)= (i-±-}g( x ) + l-U(At,x) in [0,At]xR N . 

where V e (t, x) = (1 - 6)V{t - At, x) + 9V(t, x) and t e = t - (1 - 0)At. 

From the previous results for U the existence, uniqueness, and properties of V 
easily follow. 

Proposition 9.4. Assume that (Al), (|Yip. (|Y2|) . and t/ie CFL conditions (143)) 
ZioZrf, and t/iat fc G (0, 1) and At < (2fc A 2fc 1 )~ 1 . 

(a) There exists a unique solution V € Cb{Qr) of (|9.2[) - (|9.3p . 

f&) There is a constant Ct > independent of k, At, Ax such that 

(i) |V| <Ct, 

(ii) |V(f,a;) - V(t,S/)| < C T \x-y\ for all t e [0,T], x,y,£R N , 

(in) \V{ Si ,x)-V(s2,x)\<Ct\si-s 2 \ 1/2 for all si, s 2 € [0, T], x, € M w . 

(cj Let y € Cb{Qr) o,nd V € Cb(Qr) be sub- and supersolutions of (|9.2[) — (|9.3[) 
corresponding to coefficients (o~ a ,b a , c a , f a , g) and (a a ,b a ,c a , f a ,g) respectively. 
Then there is a constant Ct > independent ofk, At, Ax such that for all t € [0, T], 

\V(t, •) - V(t, -)lo < Cr(|s - ,g|o +isu P [(|(7| A |L>| )|c Q - c a \ + \f a ~ Ho] 

+t 1 / 2 sup[\a a - a Q | + \V - b a \ }) . 

Proof. First note that the initial data on [0, At] is uniformly bounded and Lipschitz 
continuous in x and t by construction and Corollary [931 

(a) Existence of a bounded x-continuous solution follows from repeated use of 
Lemma 19.11 since we have initial conditions on [0, At]. Continuity in time follows 
from Theorem 19.21 (with x = y) since the data is t-continuous. 

(b) Part (i) and (ii) follow from Corollary 19.31 since the initial data is uniformly 
bounded and x-Lipschitz in [0, At]. To prove part (hi) we assume Si < s 2 and let 
U(t,x) and U(t,x) solve (|9.2p with data 

(o- a (t+si, x), b a {t+si,x),c a {t+si,x), f a {t+si,x), V(si,x)) and (0, 0, 0, 0, V(s u x)) 

respectively. Note that for t G [0,T— si], U(t,x) = V(si,x) and U(t,x) = V(t + 
si, x) where V is the unique solution of (|9.2[) - (|9.3p . By part (c) we then get 

|V(t + ai, •) - .)|o = ]17(*, •) - ^(*, -)lo 

<C T (0 + tsup[|r| o + |^|o|c Q |o]+t 1/2 sup[| ( 7 Q |o + |fo Q |o]) for t > 0, 

and hence part (iii) follows. 

(c) Note that by construction of the initial data and Theorem 19.21 with x = y, 
the result holds for t € [0, At], and then the result holds for any t > At by another 
application of Theorem 19.21 with x = y. □ 



20 



DEBRABANT AND JAKOBSEN 



Using Krylov's method of shaking the coefficients, we will now find smooth 
subsolutions of (|9.2|) . First we introduce the auxiliary equation 



(9.4) S At V%t,x) = inf {L%[T„ e V £ > e {t,-)](r + s,x + e) 



0<s<e 
|e|<£ 



+ c a (r + s.x + e)V e ' e {t,x) + f a (r + s,x + e)\ in (At, T] x 

J r=t e -At-e 2 

(9.5) V s (t,x) = (l-±.^g(x) + ± i V s (At,x) in [0,At] 



where r e <f)(t,x) = <j)(t,x + e) and U e (At, a;) is obtained by first solving l|9.4[) for 
discrete times f n = nAt. For this equation to be well-defined for t G (At, T], the 
data and J/Jl'^ must be defined for t S (—At — e 2 , T + e 2 ]. But this is ok since one 
can easily extend these functions to t G [— r, T + r] for any r > in such a way that 
(Al), |[YI]) . (|Y2|) still hold. Also note that 



M 



(9.6) L£[r_ e V e ' 9 (V)](r + » >a; + e) = E {v E ' e (t,x + y%+(r + s, x + e)) 

i=l 

-2^ e (t, z) + V e ' fl (t,a; + ^'"(r + s,x + e))}, 

and hence (|9.4|) is an equation of the same type as (|9.2p (with different A and 
shifted coefficients) satisfying (Al), (|YTj) , (|Y2| whenever (jjO)) does. 

By Proposition [O] there is a unique solution V £ of l[9~l |) -([9~5 ]) in [0, T + At + 

e 2 } x R N . Let [/ e (t, x) := V e {t + At + e 2 , x) and define by convolution, 

(9.7) U e (t,x)=f [ U E (t~ s,x-e)p e (s,e)dsde, 

Jm N Jo 

where e > 0, p E (t,x) = p^/o(^f), and 

pG(7 oc (R A+i ); p > Q; supp/9C [0,1] x {|.t| < 1}, /p=l. 

Note that U e is well defined on the time interval [—At, T]. By the next result it is 
the sought after smooth subsolution of (|9.2j) . 



Proposition 9.5. Under the assumptions of Proposition \97^\ the function U e de- 
fined in (|9.7p satisfies 

(i) U e e C*°°((-At,r) x R N ), \U s \i < C, \D m d?U e \ < Ce 1 -™- 271 for n,m e N. 
(«J 7/y is t/ie solution of (f9~2|) - ([93]) . t/ien |J7 e - V| < C(e + At 1 / 2 ) m Q T . 
(Hi) U s is a subsolution of (|9.2I) m Q^. 

Proof. The regularity estimates in (i) are immediate from properties of convolutions 
and the regularity of V s . The bound on U e — V (in [0, T]) in (ii) follows from 
Proposition ^. 41 (c) and (Al) which imply 

\V e -F|o < C(e + At 1 / 2 ), 

and regularity of V £ along with properties of convolutions, 

\U e - V s \ < \U e - U E \ + \V e {- + At + e 2 , ■) - V E \ < \V e \i(e + At 1 ' 2 ). 
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To see that U £ is a subsolution of (19. 2j) . first note that from the definition of U £ 
and (|9.4p it follows that 

S At U £ {t,x) < L%[T- e U £ > e (t, +s,x + e) 

+ c a (t e + s,x + e)U £ - e {t, x) + f a (t e + s,x + e) 

for all (t,x) G [-e 2 ,T] x R N , |e|,s 2 < s, and a S A Now we change variables from 
(t + s, x + e) to (i, x) to find that 

S At U e (t- «,i -e) < L£[T_ e C7 £ < e (i- S ,-)](i e ,:E) 
+c Q (i e , cc) (f - s, a; - e) + / Q (i e , x) 

for all (i, a;) € [0, T] x l w , |e|, s 2 < e, and a <E A. Then we multiply by p £ (s, e) and 
integrate w.r.t. (s,e). To see what the result is, note that 

M 
i=l 

-2U £ (t -s,x-e) + U £ (t-s,x + vt]~{r, x) - e)\, 

and hence 

J J L^ e U £ (t~ Sl -)}(r,x)p £ ( Sl e)d S de = L a [U £ (t,-)}(r,x). 

For the whole equation we then have, 

6 At U £ (t,x)< Lt [U 6 e (t,-)](t e ,x)+c a (t e , x) U 9 £ (t,x) + f a (t e , x) 

for all (t, x) € Qt and a E A. Since this inequality holds for all a, it follows that 
U £ is a subsolution of ()9.2|) in all of Qt- Q 



We are now in a position to prove the error estimate given in Theorem 

Proof of TheoremlKM Let U £ be defined in (|9T7)l . By Proposition 13] (i) and 
Lemma ETT1 (a), 

d t U e - inf (L^f/^VWtV) + c a (t e ,x)U £ > e (t,x) + f a (t e ,x)\ 

< \l^i\ d 2 Us \ oM + c{(\d 2 U £ \ + \dfU £ \ Q + \d 2 DU £ \ Q + \d 2 D 2 U £ \ Q )At 2 
+ (\DU £ \ Q + --- + \D i U £ \ )k 2 } 

< c{|l - 26»| £ - 3 A< + s~ 5 At 2 + e^ 3 fc 2 | 
in Qt- Moreover, by Proposition 19. 51 fii), 

g{x) = U(0, x) > U £ (0, x) - C(e + At 1 ' 2 ). 
It follows that there is a constant C > such that 

U £ - Ce sup ° l cQ l°*{e + At 1/2 +t(\l - 29^ At + e~ 5 At 2 + e~ 3 fc 2 ) } 
is a classical subsolution of (|1.1|) - (|1.2|) . By the comparison principle 
U £ -Ce sup «\ cC "\ ot {e + At 1/2 +t(\l-28\e- 3 At + e- 5 At 2 +e- 3 k 2 ^}<u in Q T , 
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We minimize w. r. t. 



U < 



if 
if 



Qi 



and hence by Proposition 19.51 (ii) , 

U-u = (U - U e ) + (U e - u) < c{e + Ai 1 / 2 + |l-26»|£- 3 At + £- 5 At 2 + e- 3 fc 2 }. 

and find that 

f CiAt 1 /* + k 1 / 2 ) 
\C{At x l* + k x / 2 ) 

The lower bound onw-f follows with symmetric - but much easier - arguments 
where a smooth supersolution of the equation (jl.ip is constructed. Consistency and 
comparison for the scheme (|9.2j) is then used to conclude. In view of Lemma |4.1| 
the lower bound is a direct consequence of Theorem 3.1 (a) in [2]. □ 

10. Numerical results 

In the following, we apply the LISL and MCSL schemes to linear and convex test 
problems in two space-dimensions. Hence all problems in this section are indepen- 
dent of (3. For the LISL scheme, we choose k = \J Ax and a regular triangular grid, 
whereas for the MCSL scheme we choose k = Ax and a regular rectangular grid. If 
not stated otherwise, we use 9 = (explicit methods), CFL condition At = k 2 , and 
approximation 15 . 1 1 5 1 for L Q,/3 . As error measure we will always use the L°°-norm, 



and the error rates are calculated as t; 



In ||ej 1 1 — In ||ej_i | 



All calculations are 



In || AaJiH— In \\Axt 

done in Matlab, on an INTEL Core2 Duo Mobile T7700, 2.4Ghz Laptop. 

10.1. Linear problem with smooth solution. Our first problem is taken from 
[4] and has exact solution u(t, x) = (2 — t) smxi sinx2. The coefficients in (jl.lj) are 
given by 

f a (t, x) = sin x± sin x 2 [(1 + 2/3 2 ) (2 - t) - 1] 

— 2(2 — t) cos x\ cos x% sin(a;i + x-i) cos(xi + x 2 ), 

N 



c a (t , x) =0, b a (t , x) = 0, a a (t , x) 



^2 / sin(a:i + x 2 ) P 
\ cos(a^i + x 2 ) P 



We consider 1 = 0.1 and (3 = 0. Note that in the second case, the scheme 
considered in [J] is not consistent. Table [5] gives the (spatial) errors and rates 
obtained at t = 1 applying the LISL and the MCSL scheme. As expected for 



1 = 0.1 



(3 = 



Ax 


LISL MCSL 


LISL MCSL 


error rate 


error rate 


error rate 


error rate 


3.93e-2 
1.96e-2 
9.82e-3 
4.91e-3 
2.45e-3 


3.79e-2 0.86 
1.93e-2 0.97 
9.45e-3 1.03 
4.50e-3 1.07 
2.43e-3 0.89 


1.03e-3 
2.57e-4 2.00 
6.42e-5 2.00 
1.61e-5 2.00 
4.01e-6 2.00 


3.94e-2 0.87 
1.98e-2 0.99 
9.94e-3 0.99 
4.70e-3 1.08 
2.45e-3 0.94 


1.03e-3 
2.57e-4 2.00 
6.43e-5 2.00 
1.61e-5 2.00 
4.02e-6 2.00 



Table 2: Results for the smooth linear problem at t 
grid adapted to monotonicity 



1 with (3 2 = 0.1 and (3 = 0, 



smooth solutions, in both cases we obtain order one for the LISL scheme and order 
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two for the MCSL scheme. Here, we have chosen the grid points such that the 
solution is monotone in between. Without this, we would still obtain order one for 
the LISL scheme, but no convergence for the MCSL scheme. The reason is that for 
non-monotone data, the interpolation error of monotone cubic interpolation reduces 
to second order, and so the choice k = Ax is not longer appropriate. 

10.2. Linear problem with non-smooth solution. The second problem we test 
is a problem with non-smooth exact solution in [— 7T, it] 2 given by 

x 2 fsinf- for - ir < x x < 0, 
u(t,x) = (l + t)sm—l . 2 

2 I sin ^ for < x\ < n. 

The coefficients in (jl.ip are given by 



f a (t,x) 



sin 



x 2 J sin ^- (l + ^(sin 2 x\ + sin 2 a; 2 )) for - n < xi < 
2 1 sin^- (l + ^t«(sin 2 xi +4sin 2 a;2)) for < x% < tt 



smii sm X2 cos 



x 2 J ^ cos ^- for - 7r < xi < 
T 1 ±±* cos f> for < sci < 7T 



c Q (i, a?) =0, b a (t,x) = 0, <i a (t, x) = y/2 ( Sm Xl 

and we pose Dirichlet boundary conditions. This is a monotone non-smooth prob- 
lem, and we obtain order one half applying the LISL scheme and order one applying 
the MCSL scheme, i. e. reduced rates, see Table [3] 



Aa; 


LISL 


MCSL 


error 


rate 


error 


rate 


3.90e-2 


8.75e-3 




4.19e-3 




1.96e-2 


6.I9e-3 


0.50 


2.20e-3 


0.93 


9.80e-3 


4.38e-3 


0.50 


1.12e-3 


0.97 


4.90e-3 


3.10e-3 


0.50 


5.69e-4 


0.98 


2.45e-3 


2.19e-3 


0.50 


2.86e-4 


0.99 



Table 3: Results for the non-smooth linear problem at t = 1 



10.3. Optimal control problems with smooth solutions. 

a) We test an example from [4] with exact solution u(t, x\, x 2 ) = (| — i) sinxi sinx2- 
The corresponding coefficients and control set in (11. 1 1) are 



r 



— — t J sm xi sm X2 + I — — t 



cos 2 Xi sin 2 X2 + sin 2 x\ cos 2 X2 



— 2 sin(a;i + x?) cos(a;i + x%) coscci coscc2 



c a = 0, b a = a, a a = V2 f 8 ^ 1 +X2 \), A - {a e R 2 : a\ + a\ = 1}. 

V COS(xi + X2)! 
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As er Q does not depend on a but b a does, we choose approximation 15.1141 for 
£a,p anc j thus need only about half of the number of interpolations we would 
need if we had chosen approximation 15.1151 

b) The next test problem has exact solution u(t, x\, X2) = (2 — t) sinxi sina^ and 
coefficients and control set given by 

f a (t, x) = (1 — t) sinxi sin £2 — 2aitt2(2 — t) cos x\ cos X2, 

c a {t,x)=0, b a (t,x,a) = 0, £T Q = v^f^Y A= {a eW 2 : a\+ a% = 1}. 

In both examples, the control is discretized on the unit circle by grid points. 
The results at t = 0.5 are given in Table 2] , where again the grid is adapted to 
monotonicity. As expected for smooth solutions, the LISL scheme yields a numerical 
order of convergence of one, whereas the MCSL scheme yields order two. 



Ax 


LISL MCSL 


LISL MCSL 


error rate 


error rate 


error rate 


error rate 


3.93e-2 
1.96e-2 
9.82e-3 
4.91e-3 
2.45e-3 


3.01e-2 
1.61e-2 0.91 
8.03e-3 1.00 
3.94e-3 1.03 
2.03e-3 0.96 


8.40e-4 
2.12e-4 1.98 
5.30e-5 2.00 
1.33e-5 2.00 
3.32e-6 2.00 


2.18e-2 
1.07e-2 1.03 
5.45e-3 0.97 
2.55e-3 1.10 
1.34e-3 0.92 


5.14e-4 
1.29e-4 2.00 
3.21e-5 2.00 
8.03e-6 2.00 
2.01e-6 2.00 



Table 4: Results for optimal control problems at t = 0.5, grid adapted to mono- 
tonicity 



10.4. Convergence test for a super-replication problem. We consider a test 
problem from [5J which was used to test convergence rates for numerical approxi- 
mations of a super-replication problem from finance we will consider in Subsection 
110.51 The corresponding PDE is 
(10.1) 

inf \a\ut{t,x) - ^tr (<T a (t,x)<r aT (t,x)Du(t,x)) 1 = f(t,x), 0<xi,x 2 <3 
af+a|=l L 2 J 

withcr Q (t,a;) = (^2^{^f) axsdr >( x ) = x (3~ x )- We take = l+^-e"^"^ 

as exact solution as in [BJ, and then / is forced to be 

/(*, x) = i \ u t - ~xfx 2 u XlXl - ^2(3 - x 2 ) 2 u X2X2 

-\j ^-Ut + ^xjx 2 u XlXl - ^x|(3 - x 2 ) 2 u X2X ^j + (x 1 ^/x^ 3 (3 - x 2 )u XlX2 S j 

In [BJ rj{x) — x, while we take rj(x) = x(3 — x) to prevent the LISL scheme from 
overstepping the boundaries. Note that changing r\ does not change the solutions 
as long as 77 > in the interior of the domain, see [BJ, and hence the above equation 
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is equivalent to the equation used in [6]. The initial values and Dirichlet boundary 
values at xi = and X2 — are taken from the exact solution. As in [6], at 
x = 3 and y = 3 homogeneous Neumann boundary conditions are implemented. 
To approximate the values of a\,a2, the Howard algorithm is used (see [6]), which 
requires an implicit time discretization, so we choose 9 = 1. The minimization is 
done over a\ t k + iot-2,k = e 27rlk / 2N&!C , k = 1, . .., A^az, where N^x is the number of 
space grid points, i. e. Na x = 3/ Ax. 

The results at t = 1 are given in Table [5] Again, the numerical order of conver- 
gence is approximately one when the LISL scheme is used and approximately two 
for the MCSL scheme. Note that compared to the results in 4J, for comparable 
accuracies the LISL scheme is about ten times faster, the MCSL scheme about 100 
to 1000 times faster. 



(a) LISL (b) MCSL 



Ax 


error 


rate 


time in s 


Ax 


error 


rate 


time in s 


1.50e-l 


2.01e-l 




0.71 


3.00e-l 


8.21e-2 




1.17 


7.50e-2 


9.49e-2 


1.08 


6.76 


1.50e-l 


1.83e-2 


2.16 


11.58 


3.75e-2 


4.29e-2 


1.15 


75.73 


7.50e-2 


5.03e-3 


1.86 


149.24 


1.87e-2 


1.94e-2 


1.15 


1115.39 











Table 5: Results for the convergence test for the super-replication problem at t = 1 



10.5. A super-replication problem. We apply our method to solve a problem 
from finance, the super-replication problem under gamma constraints considered in 
[6j. It consists of solving equation (|10.1j) with / = 0, Neumann boundary conditions 
and o~ a as in Subsection 110.41 and initial and Dirichlet conditions given by 

u(t, x) = max(0, 1 — xi), t = or x\ = or X2 = 0. 

The solution obtained with the LISL scheme is given in Figure [1] and coincides with 
the solution found in [Hj . It gives the price of a put option of strike and maturity 
1, and X\ and xi are respectively the price of the underlying and the price of the 
forward variance swap on the underlying. 

Appendix A. Monotonicity of solutions of (|1.1|) 

We will discuss a condition ensuring that the solution of (| 1 . 1 11 — (| 1 . 2 j) is monotone 
along some unit direction e € M. . 

(A2) Let e e R N , \e\ = 1. For all x E R N , aeA,(3eB,h>0 

a^ p (t, x + he) = a a <P(t, x), b a ^{t, x + he) = b a '^(t, x), 
c a ^{t, x + he)> c Q ^(t, x), f a ' p {t, x + he)> f a ' P (t, x), 
g(x + he) > g(x). 

Lemma A.l. Assume (Al) and (A2). If u is a viscosity solution of ljl.l|) - (|1.2|l . 

then 

u(t, x + he) ~ u(t, x) > for all h > 0, (t, x) G Qt- 
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Figure 1: Numerical solution of super-replication problem at t = 1 



Proof. Assume that u > and c a ^ < and let v(t, x) — u(t, x + he). Since v(t, x) 
satisfies (|1.1[1 - (|1.2[) at the point it, x + he), an application of (A2) shows that it 
also is a supersolution at the point (t,x). By the comparison principle u < v and 
the theorem is proved. In the general case consider w — e~ sup ° ' j ' c '°*(m + |u|o)j 
and note that w > and the corresponding c^-coefficient c"'' 3 — sup Q a |c 0, ^|o is 
non-positive. The first result then applies to w, and hence the theorem holds for 
u. □ 

Remark A.l. This result is not so far from optimal when N > 1 and the solution u 
is non-smooth (e.g. only Lipschitz continuous). To see that, we consider the linear 
case where v = u e := Du ■ e satisfies 

Vt = tr[aD 2 v] + bDv + cv + ti[a e D 2 u] + b e Du + c e u + f e in Qt 

f{t,x) 

with e-directional derivatives a e , b e , c e and f e . If (Al) holds we can conclude from 
the comparison principle that 

it e = t;>0 if />0 and u e (0,a;)>0. 

If u is non-smooth, then / is well-defined only if a e = = b e , and the condition 
that / > is essentially equivalent to assumption (A2). Of course, it is possible to 
relax (A2) if N = 1 or solutions are more regular. 

Remark A. 2. It is important to notice that the result of Lemma [A. II also holds for 
all PDEs that satisfy (A2) after (monotone) coordinate transformations. In finance 
there are many such equations, e. g. the Black-Scholes equation for a European 
option based on two stocks, 

u t = \<r\x 2 u xx + paia 2 xyu xy + ^oly 2 u yy + r{xu x + yu y ) ~ ru, t,x,y > 0, 

u(0, x, y) = max(0, K — (x + y)), x, y > 0. 

After the change of variables (x, y) = (In x, In y), this equation reduces to a constant 
coefficient equation. Since the initial data is decreasing in x and y, the same is true 
for the solution u by Lemma [All Going back to (x, y) variables, we then find that 
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u is decreasing also in x and y. (Strictly speaking we must extend u(t, •, •) to M 2 in 
a suitable way to apply Lemma lA.lj) . 

Appendix B. The proof of Theorem 19.21 



We will prove the result when k = 0. The general case can be reduced to this 
case in a standard way by considering U/Rk and U/Rk instead of U and U. We 
use doubling of variables techniques similar to those used to prove this type of 
results for equation We take 

m = \(U(0, -)-U(0, -)) + |o, 



sup 



M 2 — 4 sup [W a -a a \ 2 Q 

a 

where <j)^ denotes the positive part of q 



i(r-n + io + (it/ioA|[/| )ic« 
- \b a - b° 

and define 



W{t,x,y) 
<t>{t, V) = m a 



U(t,x)~U(t,y), 
1 



tm - 



2s 



K T tM z 



1 



+ -R kl (t)(L + tL){e + -\x- y\ 2 ) + S(\x\ 2 + \y\ 2 ), 
2 e 

^{t, x, y) = W(t, x, y) - (f>(t, x, y) - rj(l + t), 



rh= sup ip(t,x,y) = ip(t,x,y), 

x,yeR N 



for £, S, r\ > and a maximum point (t,x,y). A maximum point exists because of 
the (5-terms in <fi. We will prove that for any sequence r)i —> 0, there is another 
sequence Si ^ such that ip{ti,xi,yi) < o(l) as I — > oo. This implies Theorem [ 
when fco = 0. To see this, fix t > 0, x, y and note that for any e > 0, 



±-K T tM 2 - \R kl (t)(L Q + tL){£ + -\x - y\ 2 ) 

2£ 2 £ 



as / 



U(t, x) — U(t, y) — mo — tm — 
< ^{t u xum) + S t {\x\ 2 + \y\ 2 ) + m (l+t)< o(l) 
In this inequality we send I — > oo and choose 

e=\x-y\ Vt 1,2 M 

to find that 

U(t,x) - U(t,y) < too +tm + t 1/2 K T M + R kl (t){L +tL)\x - y\, 
and hence Theorem 19.21 follows since t > 0,x,y were arbitrary. We will not be 



explicit about the form of the 5-terms below. Their role is only to guarantee that 
the maximum is attained at a (finite) point (t, x,y), and their contribution will 
always be o(l) as S — > (see also Section 3 in pQ). 

It is enough to prove that for every r) > 0, ip(t, x, y) < o(l) as 5 — > 0. We proceed 
by contradiction assuming there is an 77 > such that lim^o i>(ti y) > 0- By the 
definition of tp we now have W(t, x,y) > and t > for all S > small enough. 
The last statement is true since 



V>(0, i, y) < toq + L Q \x - y\ - m 



Lq . 1 . _ 



y| 2 ) - 77 < 0. 



28 



DEBRABANT AND JAKOBSEN 



The rest of the proof will aim at getting a contradiction for the case t > 0. Even if we 
do not write it like that, what we show below is that ^^ t - x ^~^~^ t - x ^ < (i) — r\ 
as 5 — > 0, and this is impossible since (t, x, y) is a maximum point of ip. 
We proceed by defining the operator IP, 

M 

IP m, , -)](r, x, y) = x + y%+(r, x), y + y^(r, y)) 

»=i 

-2<j>(t, x, y) + ef,(t, x + y%7 (r, x),y + y%r (r, y)) } . 

By the definition of L£ and L£ and (|4.ip , it follows that 

U a [W(t, ;-)](r,x,y) = 2k 2 [l% [U(t, •)] (r, x) - L a k [U(t, •)] (r, y) } . 

We set A := ^ and subtract the inequalities defining U and U (see (|8.1[) and (|9.ip ) 
to find that 

+ sup {±IL a [W e (t, ., •)](*", y) + At c Q (t e , or, y)} 

+ Ati|x - ?/| + Aim for (f, x), (t, y) S Qt, 

where W^fax.y) = (1 - 0)W(* - At,a;,y) + 0W(t,x,y) and £ e = t - (1 - 0)At. 
Note that this new "scheme" is still monotone by the definition of 11" and the CFL 
condition. Hence we may replace W in the above inequality by any bigger function 
coinciding with W at (t, x, y). By the definition of m, 

W <4> + 77(1 + t) + m in AtN x R N x R N , 

and equality holds at (t, x,y). Therefore we find that 

(*) cP(t,x,y)+r)(l + t)<<j>(i-At,x,y)+r)(l + t-At) 

+ sup^U a [^ e (t,-,-)](t ,S:,y) + AtL\x~y\ + Atm. 

Here we also used the fact that n a [r;(l + 1) + rh] = and c a < 0. Moreover we can 
Taylor expand to see that 

M 

IV a m, -)](r,x,y) = £ {(Y+ +Yr) . D x ^ + (y+ +y r ) . D v 
»=i 

+ ~tr[D^ • (^ + T + YfYr T )] + hr[D 2 yy q> ■ (Y+Y+ T + YfYf T )] 
+ itr[I^,0 • (Y+Y+ T + Y+Y+ T + Yryr T + Yf Y r T )]}, 



SL SCHEMES 



2!) 



where F± = y k f(r,x) and Yf = y k f{r,y). Note that YY T = Y <g>F for Y £ R N . 
Now we use (IY2[) along with the definition of 4>, to see that 



IF[<K* )V )](r,z,y) < ii? fcl (t)(L + ^)^2fe 2 (//- !,-..,■. - //'(/■. //)!(.»•- //) 



(a a (r, x) - a a (r, y))(a a (r, x) - a a (r, y)) T 
(b a (r, x) - b a (r, y))(b a (r, x) - b a (r, y)) T ] } + o(l), 



+ fc 2 tr 
+ fc 4 tr 

as 5 — > 0. These considerations lead to the following simplification of (J*J), 

<P(t,x,y) - <P{t - At,x,y) 



V 



At 



< e-R kl (t)(L + IL){\m 2 + \k x \x - y\ 2 ) 

£ 2 2 

+ (1 - 9)-R kl (t- At)(L + (t-At)L)(hl 2 + h x \x ~ y\ 2 ) 
£ 2 2 

+ L\x — y\ + m + o(l) 

< -R kl (t)(L + tL)(-M 2 + h^x - y\ 2 ) + L\x -y\+m + o(l) := RHS, 

as S — > 0. Now we proceed to calculate 

x ,/. \ <t>{i,x,y) - <j){t - At,x,y) 
0Atq>{t,x,y) := — . 

To do that we note that 

8At(uv) = (Saiu)v + u5 A tv - At(SAtu){6Atv). 
Since SAtR^it) = fci-R/s^i) we then see that 

6 At [R kl (t) {L + tL)] = hR kl (t) (L + tL) + R kl (t)L - AtLkxR kl (t) , 
and hence 

8At4>{l x,y)=m+ \-K T -M 2 + ii2 fcl (t) [h(L + tL) +L- AthL] (e +-\x- y\ 2 ). 

lei e 

All of this leads to 

r/ < RHS - S At (f>(t,x,y) < o(l) as 5 0. 

The last inequality follows from the bound on Kt- We have our contradiction and 
the proof is complete. 
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